function savings_EWM_onedim = generate_results_capshare_onedimension(bsimax,sample,share,wave,tcost,covariate,winter_prevuse,SMC,savedir,fbase)
% This function calculates savings from one-dimensional rule under a capped 
% budget based on average baseline consumption only. The average baseline 
% consumption measure is discretized and searched over each level. 

% Inputs: 
% (1) bsimax: number of bootstrap runs
% (2) sample: cross-sectional dataset with electricity consumption, 
% covariates, and propensity scores
% (3) share: share of households allowed to treat under the budget cap
% (4) wave: indicates whether the output is wave-specific (=3,6,or 7) or
% for the pooled sample (=0)
% (5) tcost: private marginal cost of implementing the program per household
%  >0 represents cost savings; =0 represents kWh reduction 
% (6) covariate: indicates covariate for analysis: income, size, vintage,
% minimum of baseline consumption, maximum of baseline consumption, or
% standard deviation of consumption
% (7) winter_prevuse: indicates whether baseline consumption is calculated 
% as the mean of consumption in winter months (Jan and Feb) or as the mean
% of specified pre-treatment periods.
% (8) SMC: =1 use social marginal cost; =0 use retail electricity price
% (9) savedir: directory to save the savings table 
% (10) fbase: string used to indicate the propensity score and baseline
% months specification for output table 

% Output: 
% (1) savings_EWM_onedim: a table summarizing the savings point estimate,
% share of households to treat

%% Gather inputs

[X,Y,D,n,ps]=generate_inputs_quadrant_rule(sample,wave,tcost,covariate,SMC);

%% Estimate g using IPW
g = Y.*(D./ps - (1-D)./(1-ps)); % elements of W(G) - W(G_0)

%% Level search

% sort rows of X -> [Xr, gr]
[Xr, Ind] = sortrows(X);
gr = g(Ind);

% descretize baseline usage
dx1=10;
[N,boo1,xx1]=histcounts(Xr(:,2),...
    'BinWidth',[dx1],'Normalization','count'); %

nu = length(boo1);

Xu=nan(nu,1); gu=nan(nu,1); nw = nan(nu,1);
k=1;
for i = 1:length(boo1)
       bin_boo = logical((xx1==i));
       % identify individuals in each level of baseline consumption
       Xu(k,:)=boo1(i);
       gu(k,1) = sum(gr(bin_boo)); % net savings for each unique covariate combo
       nw(k,1) = sum(bin_boo); % count number of individuals in the pair of covariate combo
       k=k+1;
end

x1 = Xu(:,1);

% Create grids of midpoints between distinct values of x1/x2
u_x1 = boo1';

eps_x1 = min(u_x1(2:end) - u_x1(1:end-1))/2;
grid_x1 = [u_x1(1)-eps_x1; (u_x1(2:end) + u_x1(1:end-1))/2; u_x1(end)+eps_x1];
k1 = size(grid_x1,1);

% Grid search 
tic1= tic;
m_W = Inf*ones(k1);
v_x1 = grid_x1;
share_candidate=nan(k1,1);
fractn=nan(k1,1);

minW = Inf;

for i1=1:k1
    boo_i1=(x1 - grid_x1(i1));
        for sign1 = -1:2:1 %for i = [-1,1]
                select = ((sign1.*boo_i1)>=0); 
          
                %must meet both threshold conditions to be selected as treated
                treat_percent_candidate = sum(nw.*select)/n;
                fraction = min(1,share./treat_percent_candidate); % rationing term
                W_rationed = sum(gu.*select).*fraction; % total savings under rationing
                %
                if W_rationed<m_W(i1)
                   m_W(i1) = W_rationed; 
                   share_candidate (i1,1) = treat_percent_candidate;
                   fractn (i1,1) = fraction;
                end
                %
                if (W_rationed < minW)
                    min_i1 = i1;
                    min_sign1 = sign1;
                    minW = W_rationed;        
                end
        end
end

fprintf('Serach takes %.2f sec\n',toc(tic1))

%% Save coefficients and other variables for graphs

in_Ghat = (min_sign1*(x1 - grid_x1(min_i1))>=0);
select_fraction=fractn(min_i1,1);
percent=sum(nw.*in_Ghat)/n*select_fraction;
savings=sum(gu.*in_Ghat)/n*select_fraction;

 %% Export rule parameters for graphs 

switch winter_prevuse 
    case 1
        s_winter = '_winter';
    case 0
        s_winter = '';
end
if (tcost)>0 % cost_savings
    if (SMC)
        s_cost = 'SMC';
    else
        s_cost = 'PMC';
    end
elseif (tcost)==0
    s_cost='kwh';
end
filename_coefs=sprintf('coef_capshare_onedim_%s_%s%s_wave%1.0f.mat',covariate,s_cost,s_winter,wave);
filename = sprintf('%s_%s',fbase,filename_coefs);
fpathf = fullfile(savedir,filename);
fprintf('Saving "%s" to "%s" ... ',filename,savedir); % must use newline symbol "\n" at the end
save(fpathf,'min_i1','min_sign1',...
    'in_Ghat','v_x1','Xu','nw','gu','n','covariate','minW',...
    'percent','savings','share_candidate','fractn','select_fraction');

%% Output for table 

savings_EWM_onedim=nan(1,4);
savings_EWM_onedim=array2table(savings_EWM_onedim,'VariableNames',{'rules','Covariates','Share\ treated','Net\ cost\ changes'}); 
format short g
savings_EWM_onedim.(1) = {'Univariate'};
savings_EWM_onedim.(2) = {covariate};
savings_EWM_onedim.(3) =round(percent*100,2);
savings_EWM_onedim.(4) =round(savings,3);

end